#ifndef _SE_MATH_
#define _SE_MATH_

#include "Core/Base.h"

#include <cmath>
#include <random>
#include <chrono>
#include <algorithm>

#include <glm/glm.hpp>
#include <glm/gtc/type_ptr.hpp>
#include <glm/gtc/matrix_transform.hpp>
#define GLM_ENABLE_EXPERIMENTAL
#include <glm/gtx/quaternion.hpp>

namespace SE
{
	// Basic math types.
    using Vector2i      = glm::ivec2;
    using Vector2f      = glm::vec2;
	using Vector3f      = glm::vec3;
	using Vector4f      = glm::vec4;
	using Matrix3f      = glm::mat3;
	using Matrix4f      = glm::mat4;
    using Quaternion    = glm::quat;

	// Math definations.
	static const Float Epsilon		= 1e-6f;
	static const Float PositiveInf	= std::numeric_limits<Float>::infinity();
	static const Float NegativeInf	=-std::numeric_limits<Float>::infinity();

	static const Float PI			= 3.14159265358979323846264338327950288f;
	static const Float OneOverPI	= 1.0f / PI;
	static const Float TwoPI		= PI * 2.0f;
	static const Float HalfPI		= PI * 0.5f;

    class Math final
    {
    public:
        static Float        abs(const Float value)
        {
            return std::fabs(value);
        }
        static Float        pow(const Float value, const Float exponent)
        {
            return std::pow(value, exponent);
        }
        static Boolean      isNan(const Float f)
        {
            return std::isnan(f);
        }
        static Float        sqr(const Float value)
        {
            return value * value;
        }
        static Float        sqrt(const Float value)
        {
            return std::sqrt(value);
        }
        static Float        invSqrt(Float value)
        {
            Float xhalf = 0.5f * value;
            I32 i = *(I32*)&value;
            i = 0x5f375a86 - (i >> 1);
            value = *(Float*)&i;
            value = value * (1.5f - xhalf * value * value);
            return value;
        }
        static Boolean      realEqual(const Float a, const Float b, const Float tolerance = std::numeric_limits<Float>::epsilon())
        {
            return std::fabs(b - a) <= tolerance;
        }
        static Float        clamp(const Float v, Float min, Float max)
        {
            return std::clamp(v, min, max);
        }

        template<typename _Ty>
        static auto max(const _Ty& val1, const _Ty& val2)
        {
            return std::max(val1, val2);
        }

        template<I32 _Size, class _Ty, glm::qualifier _Q>
        static I32 maxDimension(const glm::vec<_Size, _Ty, _Q>& v)
        {
            I32 maxDimension = 0;
            for (I32 i = 1; i < _Size; i++)
                if (v[maxDimension] < v[i])
                    maxDimension = i;
            return maxDimension;
        }

        template<I32 _Size, class _Ty, glm::qualifier _Q>
        static I32 minDimension(const glm::vec<_Size, _Ty, _Q>& v)
        {
            I32 maxDimension = 0;
            for (I32 i = 1; i < _Size; i++)
                if (v[maxDimension] > v[i])
                    maxDimension = i;
            return maxDimension;
        }

        template<I32 _Size, class _Ty, glm::qualifier _Q>
        static glm::vec<_Size, _Ty, _Q> maxEachComponent(const glm::vec<_Size, _Ty, _Q>& v1, const glm::vec<_Size, _Ty, _Q>& v2)
        {
            glm::vec<_Size, _Ty, _Q> ret;
            for (I32 i = 0; i < _Size; i++)
                ret[i] = std::max(v1[i], v2[i]);
            return ret;
        }

        template<I32 _Size, class _Ty, glm::qualifier _Q, typename... Args>
        static glm::vec<_Size, _Ty, _Q> maxEachComponent(const glm::vec<_Size, _Ty, _Q>& v1, const glm::vec<_Size, _Ty, _Q>& v2, Args... args)
        {
            glm::vec<_Size, _Ty, _Q> ret = MaxEachComponent(v1, v2);
            return MaxEachComponent(ret, args...);
        }

        template<I32 _Size, class _Ty, glm::qualifier _Q>
        static glm::vec<_Size, _Ty, _Q> minEachComponent(const glm::vec<_Size, _Ty, _Q>& v1, const glm::vec<_Size, _Ty, _Q>& v2)
        {
            glm::vec<_Size, _Ty, _Q> ret;
            for (I32 i = 0; i < _Size; i++)
                ret[i] = std::min(v1[i], v2[i]);
            return ret;
        }

        template<I32 _Size, class _Ty, glm::qualifier _Q, typename... Args>
        static glm::vec<_Size, _Ty, _Q> minEachComponent(const glm::vec<_Size, _Ty, _Q>& v1, const glm::vec<_Size, _Ty, _Q>& v2, Args... args)
        {
            glm::vec<_Size, _Ty, _Q> ret = MinEachComponent(v1, v2);
            return MinEachComponent(ret, args...);
        }

        template<I32 _Size, class _Ty, glm::qualifier _Q>
        static _Ty maxComponent(const glm::vec<_Size, _Ty, _Q>& v)
        {
            _Ty maxValue = std::abs(v[0]);
            for (I32 i = 0; i < _Size; i++)
                if (maxValue < abs(v[i]))
                    maxValue = abs(v[i]);
            return maxValue;
        }

        template<I32 _Size, class _Ty, glm::qualifier _Q>
        static _Ty minComponent(const glm::vec<_Size, _Ty, _Q>& v)
        {
            _Ty minValue = abs(v[0]);
            for (I32 i = 0; i < _Size; i++)
                if (minValue > abs(v[i]))
                    minValue = abs(v[i]);
            return minValue;
        }

        template<class _Ty, glm::qualifier _Q>
        static glm::vec<2, _Ty, _Q> linearInterpolate(const glm::vec<2, _Ty, _Q>& v1, const glm::vec<2, _Ty, _Q>& v2, const Float t)
        {
            return (1 - t) * v1 + t * v2;
        }

        template<I32 _Size, class _Ty, glm::qualifier _Q>
        static glm::vec<_Size + 1, _Ty, _Q> fill(const glm::vec<_Size, _Ty, _Q>& v, const _Ty& value)
        {
            glm::vec<_Size + 1, _Ty, _Q> ret;
            for (I32 i = 0; i < _Size; i++)
                ret[i] = v[i];
            ret[_Size] = value;
            return ret;
        }

        template<I32 _Size, class _Ty, glm::qualifier _Q>
        static glm::vec<_Size - 1, _Ty, _Q> eliminate(const glm::vec<_Size, _Ty, _Q>& v)
        {
            glm::vec<_Size - 1, _Ty, _Q> ret;
            for (I32 i = 0; i < _Size - 1; i++)
                ret[i] = v[i];
            return ret;
        }

        static Float        randomFloat(const Float min, const Float max)
        {
            static std::random_device rd;
            static std::mt19937 gen(rd());
            std::uniform_real_distribution<Float> dis(min, max);
            return dis(gen);
        }
        static I16          randomInt16(const I16 min, const I16 max)
        {
            static std::random_device rd;
            static std::mt19937 gen(rd());
            std::uniform_int_distribution<int32_t> dis(min, max);
            return dis(gen);
        }
        static I32          randomInt32(const I32 min, const I32 max)
        {
            static std::random_device rd;
            static std::mt19937 gen(rd());
            std::uniform_int_distribution<int16_t> dis(min, max);
            return dis(gen);
        }
        static Double       randomDouble(const Double min, const Double max)
        {
            static std::random_device rd;
            static std::mt19937 gen(rd());
            std::uniform_real_distribution<Double> dis(min, max);
            return dis(gen);
        }
        static Vector2f     randomVector2f(const Float min, const Float max)
        {
            return Vector2f(Math::randomFloat(min, max), Math::randomFloat(min, max));
        }
        static Vector3f     randomVector3f(const Float min, const Float max)
        {
            return Vector3f(Math::randomFloat(min, max), Math::randomFloat(min, max), Math::randomFloat(min, max));
        }
        static Vector4f     randomVector4f(const Float min, const Float max)
        {
            return Vector4f(Math::randomFloat(min, max), Math::randomFloat(min, max), Math::randomFloat(min, max), Math::randomFloat(min, max));
        }
        static Matrix3f     randomMatrix3f(const Float min, const Float max)
        {
            Matrix3f result;
            for (int32_t i = 0; i < 3; ++i)
                for (int32_t j = 0; j < 3; ++j)
                    result[i][j] = Math::randomFloat(min, max);
            return result;
        }
        static Matrix4f     randomMatrix4f(const Float min, const Float max)
        {
            Matrix4f result;
            for (int32_t i = 0; i < 4; ++i)
                for (int32_t j = 0; j < 4; ++j)
                    result[i][j] = Math::randomFloat(min, max);
            return result;
        }

        template<typename _Ty>
        static void hashCombine(U32& seed, const _Ty& v)
        {
            seed ^= std::hash<_Ty> {}(v)+0x9e3779b9 + (seed << 6) + (seed >> 2);
        }

        template<typename _Ty, typename... _Tys>
        static void hashCombine(U32& seed, const _Ty& v, _Tys... rest)
        {
            hashCombine(seed, v);
            if constexpr (sizeof...(_Tys) > 1)
                hashCombine(seed, rest...);
        }
    };
} // !namespace SE

#endif // !_SE_MATH_